function slowtest(funct,N,D,M,eps)
% slowtest(funct,N,D,M,eps)
% 	Parameters:
% 		funct:  the function used to create a time series
% 		N: 	number of elements in the time series
% 		D:	inital displacement
% 		M:	kernel size
% 		eps:	perturbation used ot make displacement a function of time
% 	Outputs:
% 		Only the plots produced by this function
% 	Purpose of this code:  
% 		This is a test designed to use the fractional delay filter to shift a time series with where 
% 		displacement is a function of time.  This produces two suppression plots one where only two terms
% 		of the taylor expansion are used and compared to the same time series shifted by taking filtering at each time step.
hold off;

mesh = 1:N;
ts = funct(mesh);
ts = lowpass(ts);

slowts = slowshift(ts,D,M,eps);

D1M = D+((N./2)).*eps;


%% shift timeseries by S1 and then by S2
F11 = FDtesteps(ts,D1M,M,@sinc);
F12 = FDtesteps(ts,D1M,M,@sincdir);
F13 = FDtesteps(ts,D1M,M,@sinc2dir);


for j = 1:N
    delta(j) = (j-1-(N./2)).*eps;
end

for j = 1:N
    tsint10(j) = F11(j);
    tsint11(j) = F11(j) + delta(j).*F12(j);
    tsint12(j) = F11(j) + delta(j).*F12(j) + (1./2).*(delta(j)^2).*F13(j);
end

for j = 1:N
    Dadj = D1M - (M-1)/2; 
    Dint = round(Dadj); 
    k = (-(M-1)/2):1:((M-1)/2); 
    s = (Dadj-Dint)-k; 
    w = getWindow('blackman3',M);
    p = sinc(s).*w;
    q = sincdir(s).*w;
    q2 = sinc2dir(s).*w;
    
    intf0 = (p);
    intf1 = (p + delta(j).*q);
    intf2 = (p + delta(j).*q + (1./2).*(delta(j)).^2.*q2);
    n0(j) = sum(intf0);    
    n1(j) = sum(intf1);    
    n2(j) = sum(intf2);
end

tsint10 = tsint10./n0;
tsint11 = tsint11./n1;
tsint12 = tsint12./n2;

tsd1 = tsint11;
tsd2 = tsint12;

%plot(mesh,slowts,'r');
%hold on;
%plot(mesh,tsd1,'g');

[psorig,freq] = mypsd(ts);
err = abs(slowts-tsd1);
errps = mypsd(err);
suppr = sqrt(psorig./errps);
figure;
loglog(freq,suppr,'r');
hold on;

err2 = abs(slowts-tsd2);
errps2 = mypsd(err2);
suppr2 = sqrt(psorig./errps2);
loglog(freq,suppr2,'g');

end

